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Abstract We present a numerical method to compute expectations of func- 
tionals of a piecewise-deterministic Markov process. We discuss time dependent 
functional as well as deterministic time horizon problems. Our approach is 
based on the quantization of an underlying discrete-time Markov chain. We 
obtain bounds for the rate of convergence of the algorithm. The approxima- 
tion we propose is easily computable and is flexible with respect to some of the 
parameters defining the problem. Two examples illustrate the paper. 
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1 Introduction 



The aim of this paper is to propose a practical numerical method to approximate 
some expectations related to a piecewise-deterministic Markov process thanks 
to the quantization of a discrete-time Markov chain naturally embedded within 
the continuous-time process. 

Piecewise-deterministic Markov processes (PDMP's) have been introduced 
by M.H.A. Davis in [5] as a general class of stochastic models. PDMP's are a 
family of Markov processes involving deterministic motion punctuated by ran- 
dom jumps. The motion depends on three local characteristics namely the 
flow the jump rate A and the transition measure Q, which specifies the post- 
jump location. Starting from the point x, the motion of the process follows the 
flow <^(x,t) until the flrst jump time Ti, which occurs either spontaneously in a 
Poisson-like fashion with rate A(<I>(a;, t)) or when the flow t) hits the bound- 
ary of the state space. In either case, the location of the process at the jump 
time Ti, is selected by the transition measure (3(<i>(a;, Ti), •) and the motion 
restarts from this new point denoted Zi. We define similarly the time S2 
until the next jump, T2 = Ti + S2 with the next post-jump location defined by 
Z2 = and so on. Thus, associated to the PDMP we have the discrete-time 
Markov chain (Z„, Sn)n&ii given by the post-jump locations and the inter-jump 
times. A suitable choice of the state space and the local characteristics A and 
Q provides stochastic models covering a great number of problems of operations 
research as described in [5] section 33. 

We are interested in the approximation of expectations of the form 



where {Xtj^^^ is a PDMP and I and c are some non negative, real-valued, 
bounded functions and dE is the boundary of the domain. Such expectations 
are discussed by M.H.A. Davis in [5], chapter 3. They often appear as "cost" 
or "reward" functions in optimization problems. The first term is referred to 
as the running cost while the second may be called the boundary jump cost. 
Besides, they are quite general since M.H.A. Davis shows how a "wide variety of 
apparently different functionals" can be obtained from the above specific form. 
For example, this wide variety includes quantities such as a mean exit time and 
even, for any fixed t > 0, the distribution of Xt (i.e. E2,[li?(X4)] where is a 
measurable set). 

There are surprisingly few works in the literature devoted to the actual com- 
putation of such expectations, using other means than direct Monte Carlo simu- 
lations. M.H.A Davis showed that these expectations satisfy integro-differential 
equations. However, the set of partial differential equations that is obtained is 
unusual. Roughly speaking, these differential equations are basically transport 
equations with a non-constant velocity and they are coupled by the boundary 
conditions and by some integral terms involving kernels that are derived from 
the properties of the underlying stochastic process. The main difficulty comes 
from the fact that the domains on which the equations have to be solved vary 
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from one equation to another making their numerical resolution highly problem 
specific. Another similar approach has been recently investigated in [HIT]. It is 
based on a discretization of the Chapman Kolmogorov equations satisfied by the 
distribution of the process (Xt)j>g. The authors propose an approximation of 
such expectations based on finite volume methods. Unfortunately, their method 
is only valid if there are no jumps at the boundary. Our approach is completely 
different and does not rely on differential equations, but on the fact that such 
expectations can be computed by iterating an integral operator G. This op- 
erator only involves the embedded Markov chain (Z„,S'„)„gN and conditional 
expectations. It is therefore natural to propose a computational method based 
on the quantization of this Markov chain, following the same idea as [B]. 

There exists an extensive literature on quantization methods for random 
variables and processes. The interested reader may for instance consult [8j, 
[S] and the references within. Quantization methods have been developed re- 
cently in numerical probability or optimal stochastic control with applications 
in finance (see e.g. [T], and [H]). The quantization of a random variable X 
consists in finding a finite grid such that the projection X of X on this grid 
minimizes some norm of the difference X ~ X . Roughly speaking, such a 
grid will have more points in the areas of high density oi X. As explained for 
instance in [Hj, section 3, under some Lipschitz-continuity conditions, bounds 
for the rate of convergence of functionals of the quantized process towards the 
original process are available. 

In the present work, we develop a numerical method to compute expecta- 
tions of functionals of the above form where the cost functions I and c satisfy 
some Lipschitz-continuity conditions. We first recall the results presented by 
M.H.A. Davis according to whom, the above expectation may be computed by 
iterating an operator denoted G. Consequently, it appears natural to follow the 
idea developed in [6, namely to express the operator G in terms of the under- 
lying discrete-time Markov chain (Z„, Sn)nen and to replace it by its quantized 
approximation. Moreover, in order to prove the convergence of our algorithm, 
we replace the indicator function l{x _ £dE} contained within the functional by 



some Lipschitz continuous approximation. Bounds for the rate of convergence 
are then obtained. However, and this is the main contribution of this paper, we 
then tackle two important aspects that had not been investigated in [B] . 

The first aspect consists in allowing c and / to be time depending functions, 
although still Lipschitz continuous, so that we may compute expectations of the 



This important generalization has huge applicative consequences. For instance, 
it allows discounted "cost" or "reward" functions such as l{x,t) = e~^*l{x) and 
c{x,t) = e~'^*c(a;) where 5 is some interest rate. To compute the above expec- 
tation, our strategy consists in considering, as it is suggested by M.H.A. Davis 
in [5], the time augmented process Xt = {Xt,t). Therefore, a natural way to 
deal with the time depending problem is to apply our previous approximation 



form 
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scheme to the time augmented process {Xt)t>Q. However, it is far from obvious, 
that the assumptions required by our numerical method still hold for this new 
PDMP {Xt)t>o- 

The second important generalization is to consider the deterministic time 
horizon problem. Indeed, it seems crucial, regarding the applications, to be 
able to approximate 



/ l{Xt,t)dt+ c{X^-,T,)l{x 

•^0 T,<t, 



for some fixed tf > regardless of how many jumps occur before this de- 
terministic time. To compute this quantity, we start by choosing a time N 
such that P (T/v < tf) be small so that the previous expectation boils down 



to Ex 



Jo" K^t,t)l{t<tf}dt + ciXj,- , Tj)l[x edE}'^{T,<tf} 



At first 



sight, this functional seems to be of the previous form. Yet, one must recall 
that Lipschitz continuity conditions have been made concerning the cost func- 
tions so that the indicator functions l{.<f^.} prevent a direct application of the 
earlier results. We deal with the two indicator functions in two different ways. 
On the one hand, we prove that it is possible to relax the regularity condition 
on the running cost function so that our algorithm still converges in spite of the 
first indicator function. On the other hand, since the same reasoning cannot 
be applied to the indicator function within the boundary jump cost term, we 
bound it between two Lipschitz continuous functions. This provides bounds for 
the expectation of the deterministic time horizon functional. 

An important advantage of our method is that it is flexible. Indeed, as 
pointed out in [T], a quantization based method is "obstacle free" which means, 
in our case, that it produces, once and for all, a discretization of the process in- 
dependently of the functions / and c since the quantization grids merely depend 
on the dynamics of the process. They are only computed once, stored off-line 
and may therefore serve many purposes. Once they have been obtained, we are 
able to approximate very easily and quickly any of the expectations described 
earlier. This flexibility is definitely an important advantage of our scheme over 
standard methods such as Monte-Carlo simulations since, with such methods, 
we would have to run the whole algorithm for each expectation we want to com- 
pute. This point is illustrated in Section |6] where we easily solve an optimization 
problem that would be very laboriously handled by Monte-Carlo simulations. 

The paper is organized as follows. We first recall, in Section |2] the definition 
of a PDMP and state our assumptions. In Section |3] we introduce the recursive 
method to compute the expectation. Section |4] presents the approximation 
scheme and a bound for the rate of convergence. The main contribution of the 
paper lies in Section [5] which contains the generalizations to the time dependent 
parameters and the deterministic time horizon problems. Eventually, the paper 
is illustrated by two numerical examples in Section [6] and concluded in Section [7] 
while technical results are postponed to the Appendix. 
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2 Definitions and assumptions 



For all metric space E, we denote B{E) its Borel cr-field and B{E) the set of 
real- valued, bounded and measurable functions defined on E. For a,b G M, 
denote a Ab = min(a, 6), a V 6 = max(a, b) and a"*" = a V 0. 



Definition of a PDMP 

In this first section, let us define a piecewise-deterministic Markov process and 
introduce some general assumptions. Let M be a finite set called the set of the 
modes that will represent the different regimes of evolution of the PDMP. For 
each rn G M, the process evolves in an open subset of R'^. Let 

E={{m,C),meM,CeE^}. 

This is the state space of the process {Xt)t(£R+ — {mt,Ct)teR+ ■ Let dE be 
its boundary and E its closure and for any subset Y of E, V denotes its 
complement. 

A PDMP is defined by its local characteristics {^m, Am, Q,n)meM- 

• For each m E M, <i>„i : M'^xM — > M'' is a continuous function called the fiow 
in mode m. For all i G M, <i>m(', t) is an homcomorphism and t — > ^m{', t) 
is a semi-group i.e. for all ( G M'', <i>m(C, t+ s) — ^„i{^m{C^ s),t). For all 
X — (to, C) G E, define now the deterministic exit time from E : 

t*{x) = inf{t > such that «'m(C,i) e dEm}. 

We use here and throughout the whole paper the convention inf $ = +oo. 

• For all m G Af , the jump rate A„i : E„i — > is measurable and satisfies : 

V(m, C) e i?, 3e > such that / A^ t))dt < +oo. 

Jo 

• For all m G A/, Qm is a Markov kernel on {B{E), Em) which satisfies : 

vcg:e„„ g„(c,{(m,c)r)-i. 



From these characteristics, it can be shown (see [S]) that there exists a filtered 
probability space {Q,J^,J^t,{Px)x£E) on which a process {Xt)t£R+ is defined. 
Its motion, starting from a point x £ E, may be constructed as follows. Let Ti 
be a nonnegative random variable with survival function : 

. , r e-A(-'*) ifO<t<r(x), 

where for x = (to, C) G E and t G [0,t*(a;)], 

t 



A(a;,<) = / A„($,„(C,s))ds. 
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One then chooses an i5- valued random variable Z\ according to the distribution 
Qm(*m(C,Ti), •)• The trajectory of Xt for i < Ti is : 

r (m,$™(C,i)) ift <ri, 
* \ Zi ifi^Ti. 

Starting from the point X^^ — Z\ , one then selects in a similar way S2 — T2 — Ti 
the time between Ti and the next jump, Z2 the next post-jump location and 
so on. M.H.A. Davis shows, in [5 , that the process so defined is a strong 
Markov process {Xt)tyQ with jump times (r„)„gN (with Tg = 0). The process 
(0n)neN = {Zn, Sn)ne'N whcrc Z„ = Xt„ IS the post-jump location and Sn = 
Tn — Tn-i (with 5*0 = 0) is the n-th inter-jump time is clearly a discrete-time 
Markov chain. 

The following assumption about the jump-times is standard (see for example 
[5], section 24) : 

Assumption 2.1 For all {x,t) e E x M+, E^; [^j, < +00. 

It implies in particular that Tk goes to infinity a.s. when k goes to infinity. 



Notation and assumptions 

For notational convenience, any function h defined on E will be identified with 
its component functions h„i defined on i?,„. Thus, one may write 

h{x) — hm{C) when x = (to, (^) G E. 

We also define a generalized fiow <i> : _E x M+ E such that 

^{x, t) = (to, ^>„i(C, t)) when x = (to, Q) e E. 



Define on E the following distance, for x = {m, C) and x' — (to', C) S E, 

|.-x'| = ( +°°., (1) 

[ IC — C I otherwise. ^ ^ 

For any function w in B(E), introduce the following notation 

Qw{x) ^ / w{y)Q{x,dy), = sup |w(x)|, 

and for any Lipschitz continuous function w in B{E), denote [w]^ , or if there 
is no ambiguity [w], its Lipschitz constant: 

[wf= sup 

x^yeE \x~y\ 

with the convention — = 0. 

00 

Remark 2.2 For w e B{E) and from the definition 0/ the distance on E, one 
has [w] = max„i(zM[wm]. 
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Definition 2.3 Denote Lc(-E) the set of functions w G B{E) that are Lipschitz 
continuous along the flow i.e. the real-valued, hounded, measurable functions 
defined on E and satisfying the following conditions: 

• For allx G E, w{'^{x, •)) ; [0, t*{x)[-^ R is continuous, \\m^^^*^^^~^ w{'^{x, t)) 
exists and is denoted w{^{x,t*{x))), 

• there exists [w]f e M"*" such that for all x,y G E and t e [0, t*{x) A t*{y)], 
one has: 

\wmx,t))~wmy,t))\ < [w]f\x-y\, 

• there exists [w]2 G such that for all x £ E and t,u & [0,t*{x)], one 
has: 

\w{^{x,t)) - w{^{x,u))\ < [w]^\t-u\, 

• there exists [w]f G M+ such that for all x,y G E, one has: 

\wmx,t*ix))) - wi<^iy,t*{y)))\ < [w]f\x - y\. 

Denote also Lc{dE) the set of real-valued, bounded, measurable functions defined 
on dE satisfying the following condition: 

• there exists [w]^^ G 1R+ such that for all x,y G E, one has: 

\wmx,t*{x))) - wmy,t*{y)))\ < [Hf^k - vl- 

Remark 2.4 When there is no ambiguity, we will denote instead of [w]f 
for i G {1,2,*} and [w]* instead of [w]^^ . 

Remark 2.5 In the above definition, we used the generalized flow for notational 
convenience. For instance, the definition of [w]i is equivalent to the following: 
for all m G M, there exists € such that for all G E^ and 

t G [0,t*(m,C) Ai*(m,C')]j one has: 

Let [w]i = max meM[wm\i- 

Definition 2.6 For all u > 0, denote L"(i?) the set of functions w G B{E) 
Lipschitz continuous along the flow until time u i.e. the real-valued, bounded, 
measurable functions defined on E and satisfying the following conditions: 

• For all x G E, w{^{x, ■)) : [0,t*{x)Au[-^ R is continuous and if t*{x) < u, 
then limt^t*(2.) w(<i>(x, i)) exists and is denoted w(^{x,t* (x))) , 

• there exists [w]f'" G M"*" such that for all x,y G E and t G [0,t*{x) A 
t* (y) A u] , one has: 

\wmx,t))-wmy,t))\<[w]f'^\x-y\, 

• there exists [w]2'^ G M"*" such that for all x G E and t,t' G [0,t*{x) A u], 
one has: 

\wmx,t))-wmx,t'))\<[w]^'^\t-t'\, 
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• there exists [w]f'^ E M'^ such that for all x,y € E, if t*{x) < u and 
t*{y) ^ u, one has: 

Iw^xXix))) - wiHy, t*{ym < [w]f'-\x - y\. 

Remark 2.7 For all u < u' , one has h^' (E) C L^iE) with Hf < Hf 
where i G {1, 2, *}. 



Remark 2.8 Note that Definitions \ 2.3\ and \2.6\ correspond respectively to the 
Lipschitz and local Lipschitz continuity along the flow that is, along the trajec- 
tories of the process. They can he replaced by (local) Lipschitz assumptions on 
the flow <i>, t* and w in the classical sense. 

We will require the following assumptions. 

Assumption 2.9 The jump rate A is bounded and there exists [A]i G such 
that for all x,y E E and t € [0, t*(x) A t*{y)\, one has: 

|A(<i>(x,t))-A(<f>(y,t))| < [\]i\x-y\. 

Assumption 2.10 The deterministic exit time from E, denoted t* , is assumed 
to be bounded and Lipschitz continuous on E. 

Remark 2.11 Since the deterministic exit time t* is bounded by Cf , one may 
notice that 'L'^{E) for u > Cf is no other than Lc(-B). 

Remark 2.12 Ln most practical applications, the physical properties of the sys- 
tem ensure that either t* is bounded, or the problem has a natural finite deter- 
ministic time horizon tf. Ln the latter case, there is no loss of generality in 
considering that t* is bounded by this deterministic time horizon. This leads to 
replacing Cf by tf. An example of such a situation is presented in an industrial 
example in Section \67^ 

Assumption 2.13 The Markov kernel Q is Lipschitz in the following sense: 
there exists [Q] G M+ such that for all u > and for all function w € L"(i?), 
one has: 

1. for all x,y e E and t G [0,t*{x) A t*{y) A u[, 

\Qwmx,t))-Qw{<i>{y,t))\ < [Q]Hf'"k-y|- 

2. for all x,y £ E such that t*{x) V t*{ii) < u, 

\Qw{^{x,t*{x)))-Qw{<^{y,t*{y)))\ < [Q](Hf « + - y|. 



Remark 2.14 Notice that assumption 2.13 is slightly more restrictive that its 
counterpart in JB^ (assumption 2.5) because of the introduction of the state space 
L"(i?). This is to ensure that the time augmented process still satisfies a similar 
assumption, see Section\5.1\ 
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3 Expectation 



Prom now on, we will assume that Zq = x a.s. for some x € E. For all 
fixed N GW,we intend to numerically approximate the quantity 



Jn{1,c){x) = E^; 



Jo 



(2) 



is 



where / e B{E), c e B{dE) and X^- is the left limit of Xf Thus, X^ 
the j-th pre-jump location. Since the boundary jumps occur exactly at the 
deterministic exit times from E, one has, 



Jn{1,c){x) = 



/ z(x*)dt + ^c($(z,_i,r(z,_i)))i 



{Sj=t*(Zj_i)} 



In many applications, Jn{1,c){x) appears as a "cost" or a "reward" function. 
The first term, that depends on is called the running cost and the second one, 
that depends on c, is the boundary jump cost. 

The rest of this section is dedicated to finding a formulation of the above 
expectation that will allow us to derive a numerical computation method. The 
Lipschitz continuity property will be a crucial point when it will come to proving 
the convergence of our approximation scheme. For this reason, the first step 
of our approximation is to replace the indicator function in J]\[{l,c){x) by a 
Lipschitz continuous function. Then, we will present a recursive method yielding 
the required expectation. This recursive formulation will be the basis of our 
numerical method. 



3.1 Lipschitz continuity 

We introduce a regularity assumption on I and c. 

Assumption 3.1 We assume that I G Lc(£') and c G LddE). 

Moreover, we replace the indicator function in Jm{1,c){x) by a Lipschitz 
continuous function denoted 5"^ with A> Q. Let then 



Ji^(Z,c)(a;)=E, 



/ i{Xt)dt + 5,-) 



where 5^ is a triangular approximation of the indicator function. It is defined 
on £■ X R by 

( A{t-{e{x)-\)) iortG[t*{x)-\;e{x)], 
5''{x,t) = l -A{t-{t*{x) + \)) ioxt&[t*{x);t*{x) + \\, 
y otherwise. 

For all X G E, the function 6^{x,t) goes to l^f^f,(^x^y when A goes to in- 
finity. The following proposition proves the convergence of J^{l,c){x) towards 
Jn{1,c){x) with an error bound. 
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Proposition 3.2 For all x € E , A > 0, N e N* , I £ Lc{E) and c e Lc{dE), 
one has 



\j^{l,c){x)-JNil,c)ix)\ < 



A ■ 



Proof For all a; e E, one has 

\j^{l,c){x) ~ Jn{1,c){x)\ 



E, 



N 



^c(<i>(z,_i,r(z,_i)))(<5^(z,_i,5,) ~ i{s,=t.(z,_,)}) 

r 

N 

- C'c^Ej, E \S'^{Zj^i,Sj) - l{s,=fiZj-i)}\ 



N 



We recall that the conditional law of Sj with respect to Zj-i has density s 
A($(Zj_i,s))e-^(^j-i'^) on [0;t* {Zj^i)[ and puts the weight e-^(^j-i'**(^j-i)) 
on the point t*{Zj_i). We also recall that A is bounded thanks to assump- 
tion (STOl Eventually, one has 



5^{Zj^i,s)X{^{Z,.i,s))e-^'-^'-'''Us 



\j^{l,c){x) - JNil,c){x)\ 



N 

< Cc Ea: 
J = l 
NCcCx 



t*(2,-i)-i 



< 



A 



Hence the result. 



□ 



Consequently to this proposition, we consider, from now on, the approxima- 
tion of J^{l,c){x) for some fixed A, large enough to ensure that the previous 
error is as small as required. The suitable choice of A will be discussed in 
section 14.21 



3.2 Recursive formulation 

M.H.A. Davis shows in |S], section 32, that the expectation J^{l,c){x) we are 
interested in is obtained by merely iterating an operator that we will denote G. 
The rest of this section is dedicated to presenting this method from which we 
will derive our approximation scheme in Section [4j 

Definition 3.3 Introduce the functions L, C and F defined for all x ^ E and 

t^[Q-t*{x)] by 

L{x,t) = / l[<^{x,s))ds, 
Jq 

C{x,t) = c{<i>{x,t*{x)))S'^{x,t), 
F{x,t) = L{x,t) + C{x,t), 
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along with the operator G : B{E) — >■ B{E) 



Gw;(x) = E, [F{x,Si)+w{Z^)]. 

Definition 3.4 Introduce the sequence of functions {vk)a<k<N B{E) defined 
as follows: 

vn{x) = 0, 
Vkix) = Gvk+iix), 



M.H.A. Davis then shows in [5J, equation 32.33, that for all k E {0, ...,N}, 
Thus, the quantity J^{1, c){x) we intend to approximate is none other than V(){x). 



VN-kix) = Ex 



r l{Xt)dt + ^ c(<i>(Z,_i, r 5,) 



Notice that, thanks to the Markov property of the chain (Z„,S'„)„gN, one 
has for all k G {0,...,iV- 1}, 

Gw{x) = E [FiZk,Sk+i) + wiZk+i)\Zk = x] . (3) 

Hence, for all k E {0, ...,N}, let Vk — Vk{Zk) so that one has 

Vn =0, 

Vk ^-E[F{Zk,Sk+i) + Vk+i\Zk]. 
This backward recursion provides the required quantity 

Vo = J^il,c){x). 

Consequently, we need to approximate the sequence of random variables 
{Vk)o<k<N ■ This sequence satisfies a recursion that only depends on the chain 
{Zki Sk)o<k<N ■ Therefore, it appears natural to propose an approximation 
scheme based on a discretization of this chain {Zk, Sk)Q<k<N, called quanti- 
zation, similarly to the ideas developed in [B] and [3]. 

4 Approximation scheme 

Let us now turn to the approximation scheme itself. We explained in the pre- 
vious section how the expectation we are interested in stems from the itera- 
tion of the operator G that only depends on the discrete-time Markov chain 
{Zk, Sk)a<k<N- The first step of our numerical method is therefore to discretize 
this chain in order to approximate the operator G. 

4.1 Quantization of the chain {Z/^, Sk)k<N 

Our approximation method is based on the quantization of the underlying dis- 
crete time Markov chain {Qk)k<N = {Zk, Sk)k<N- This quantization consists in 
finding an optimally designed discretization of the process to provide for each 
step k the best possible approximation of Ok by a random variable Ok which 
state space has a finite and fixed number of points. Here, optimal means that 
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the distance between Qk and 0^ in a suitably chosen LP norm is minimal. For 
details on the quantization methods, wc mainly refer to [5] but the interested 
reader can also consult [T], [2J and the references therein. 

More precisely, consider X an M'^- valued random variable such that ||^||p < 
oo and let K he a fixed integer. The optimal Lp-quantization of the random 
variable X consists in finding the best possible Lp-approximation of X by a 
random vector X taking at most K values: X £ {x^, . . . ,x^}. This procedure 
consists in the following two steps: 

1. Find a finite weig hted grid F C M« with F = {a;\ . . . , x^}. 

2. Set X = X^ where X^ = projr{X) with projr denotes the closest neigh- 
bour projection on F. 

The asymptotic properties of the ip-quantization are given by the following 
result, see e.g. [S]. 

Theorem 4.1 // E[|X|^'+''] < +oo for some 77 > then one has 

lim min \\X-X^\\p = J„ „ / \h\'^/'-'^+P\u)du, 

K^oc \r\<K P J I I 

where the law of X is Px{du) ~ h{u)\q{du) + v with v _L A^, Jp,q a constant 
and Xq the Lebesgue measure in M.'' . 

Remark that X needs to have finite moments up to the order p + 77 to ensure 
the above convergence. In this work, we used the CFVQ quantization algorithm 
described in [T], Section 3. 

There exists a similar procedure for the optimal quantization of a Markov 
chain {Xk\k£'H- There are two approaches to provide the quantized approxima- 
tion of a Markov chain. The first one, based on the quantization at each time 
k of the random variable X^ is called the marginal quantization. The second 
one that enhances the preservation of the Markov property is called Marko- 
vian quantization. Remark that for the latter, the quantized Markov process is 
not homogeneous. These two methods are described in details in [HI section 3]. 
In this work, we used the marginal quantization approach for simplicity reasons. 

The quantization algorithm provides for each time step < k < N a finite 
grid Tk of E X M"*" as well as the transition matrices (Qfc)o<fc<Ar-i from to 
Ffc+i. Let p > I such that for all k < N, and Sk have finite moments 
at least up to order p and let projr,. be the closest-neighbor projection from 
E X IR+ onto Ffe (for the distance associated to norm p) . The quantized process 
(0fc)/c<jv = {Zk, Sk)k<N takes values for each k in the finite grid F^ of £^ x M+ 
and is defined by 

{Zk,Sk) ^ projvAZk^Sk). (4) 

Moreover, we also denote respectively Ff and Ff the projections of F^: on 
E and M+. 
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Some important remarks must be made concerning the quantization. On 
the one han d, th e optimal quantization has nice convergence properties stated 
by Theorem 4.1 Indeed, the i^-quantization error \\Qk — ©fellp goes to zero 
when the number of points in the grids goes to infinity. However, on the other 
hand, the Markov property is not maintained by the algorithm and the quan- 
tized process is generally not markovian. Although the quantized process can 
be easily transformed into a Markov chain (see [HI), this chain will not be homo- 
geneous. It must be pointed out that the quantized process {&k)kei'i depends 
on the starting point Go of the process. 



In practice, we begin with the computation of the quantization grids which 
merely requires to be able to simulate the process. This step is quite time- 
consuming, especially when the number of points in the quantization grids is 
large. However, the grids are only computed once and for all and may be 
stored off-line. What is more, they only depend on the dynamics of the process, 
not on the cost functions I and c. Hence, the same grids may be used to 
compute different expectations of functionals as long as they are related to the 
same process. Our schemes are then based on the following simple idea: we 
replace the process by its quantized approximation within the operator G. The 
approximation is thus obtained in a very simple way since the quantized process 
has finite state space. 



4.2 Approximation of the expectation and rate of conver- 
gence 

We now use the quantization of the process (Ofc)fc<Ar = (^fc, Sk)k<N- In order to 
approximate the random variables (yk)k<Ni we introduce a quantized version of 
the operator G. Notice that the quantized process is no longer an homogeneous 
Markov chain so that we have different operators for each time step k. Their 
definitions naturally stem from a remark made in the previous section: recall 
that for all k e {1, N} and x £ E, 

Gw{x) = E [F{Zk-i, Sk) + w{Zk)\Zk-i = x] 



Definition 4.2 For all k G {1, N}, w e S(rf ) and z e rf_^, let 



Gkw{z)^-E F{z,Sk) + w{Zk)\Zk-i = z 



,N -1} and ze Tf. 



we then introduce the functions (vk)a<k<N ^ 

wjv(z) = 0, for all z G F^, 

Vk{z) = G'fc+iUfc+i(z), for all k e {0,.. 

Eventually, for all k G {0, N}, let 

Vk = Vk{Zk). 

Remark 4.3 The conditional expectation in Gkw{z) is a finite sum. Thus, the 
numerical computation of the sequence (V/c)fe will be easily performed as soon as 
the quantized process {Qk)k<N has been obtained. 
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Remark 4.4 We have assumed that Zq = x a.s. Thus, the quantization al- 
gorithm provides that Zq = x a.s. too. Consequently, the random variable 
Vb = vq{Zq) is, in fact, deterministic. 

The following theorem states the convergence of Vq towards Vq — J^{1, c){x) 
and provides a bound for the rate of convergence. 

Theorem 4.5 For all k e {0,...,iV}, one has v^. e hc{E). Moreover, the 
approximation error satisfies: 

\J^il, c)ix) -Vol <eNil,c,X, A) 

where 



En 



k=0 



{2[vk] + [Fl)\\Zk - ZkWp + [FUSk+i - Sk+i\\p] + 



with 

[Fh ^Ct4l]i + [cl + A[tlC,, 

[F]2 ^Q+ACc 

[vn]i <e^"^^ (^K{A,Vn-i) + nCt>[Mi{Ct'Ci + C,>j')+Ct4l]i, 

[vnh < e^"^^ (Ct^CiCx + 2Ci + CxC + {2n - l)Cx{Ct>Ci + C,)) + 

[Vn\* < [Vn]l + [t*][Vn]2- 
K] <K{A,Vn-l), 

and for all w G hdE), K{A,w) = Ei + E2A + E3[w]i + E^C^ + where 
eventually 

El = 2[/]iCt. + Ci{[t*] + 2Cl [A]i) + [c],(l + Ct^Cx) 

+Cc(2[A]iCt. +CaC2 [A]i +2[r]CA), 
E2^ CcCt'Cx[t% 
E3^ {l + Ct^Cx)[Q], 
E^^ 2C\[t*] + CtAM'^ + Ct'Cx). 



The choice of A Proposition |3 . 2| suggests that A should be as large as possi- 
ble. However, the constants [F]i, [F]2 and [u„] that appear in the bound of the 
approximation error proposed by the above theorem |4.5| grow linearly with A. 
Thus, in order to control this error, it is necessary that the order of magnitude 
of the quantization error ||Ofc — Qk\\p be at most ^. 

The convergence of the approximation scheme can be derived from theo- 
Indeed, on the one hand, one must remind that Vq = J^i^^ c)(a;) is the 



4.5 



expectation we intended to approximate and on the other hand, ||0fe — Ofc||p may 
become arbitrarily small when the number of points in the quantization grids 
goes to infinity (see e.g. [S]). An outline of the proof is presented in Appendix [C| 
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5 Time depending functionals 

We now turn to the main contribution of this paper and present two gener- 
ahzations of the previous problem. On the one hand, we will consider time 
depending functionals of the form 

Ex 

where I and c are Lipschitz continuous functions. On the other hand, we wish 
to replace the random time horizon Tn by a deterministic one denoted tf i.e. 

We will reason as follows. As it is suggested by M.H.A. Davis in [5^, we will 
introduce a transformation {Xt)t>o of the initial process {Xt)t>o by including 
the time variable into the state space i.e. (Xt) — {Xt,t). Indeed, we will see 
that both the expectation of the time depending functional and the one with 
deterministic time horizon are no other than expectations of time invariant func- 
tionals for the time augmented process {Xt)t>o- We therefore intend to apply 
the previously exposed approximation scheme to this new PDMP. However, it 
is far from obvious that the Lipschitz continuity assumptions |2.9| |2.13| and |2.10| 
still hold for this new process. 

Thus, the rest of this section is organized as follows. First, we will recall 
the precise definition of the time augmented process and prove that it satisfies 
the Lipschitz continuity assumptions required by our approximation scheme. 
Then, we will see that the time depending functional case corresponds to a 
time invariant functional for the new transformed process and may therefore 
be obtained thanks to the earlier method. Eventually, we will consider the 
deterministic time horizon problem that features an additional hurdle namely 
the presence of non Lipschitz continuous indicator functions. 

5.1 The time augmented process 

M.H.A. Davis suggests, in [5,, section 31, that the case of the time dependent 
functionals may be treated by introducing the time variable within the state 
space. Thus, it will be possible to apply our previous numerical method to the 
time augmented process. However, and this is what we discuss in this section, 
it is necessary to check whether the Lipschitz continuity assumptions still hold. 
We first recall the definition of the time augmented process given by M.H.A. 
Davis. 

Definition 5.1 Introduce the new state space 

E = ExR+ 

equipped with the norm defined by: for all ^ — (a;,t), ^' — {x' ,t') G E, let 

ie-e'i = |x-x'i + |t-i'i (5) 
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/ l{Xut)dt + y^c{X^-,T,)l^x 



l{Xt,t)dt- 



Ti<t 



c{Xj,-^,Tj)l[x^.edE} 



where the norm on E is given by ([T]) . On this state space, we define the process 

Xt^{Xt,t). 

The local characteristics of the PDMP {Xt)t>o, denoted {X,Q,^) are given for 
all ^^{x,t) eE by 



$(^,s) = ($(x,s),t + s) fors<r(a;), 
qIc,Ax {t}) ^ Q{x,A) for all AeB{E). 

Moreover, we naturally define for all £_ = {x,t) E E 

?(^) = inf{s > such that $(^, s) e dE} = t*{x) 
and 



Clearly, Assumptions 



2.9 



2.10 



still hold with [A]i = [A]i and [t*] = [t*]. 

However, proving assumption [2TT3 is more intricate. We start with the following 
lemma. 

Lemma 5.2 Let u,t > and w G L"(i?). Denote wt the function of B{E) 
defined by Wt — w{-,t). One has then Wt G L*'^"(i?) with [wtjf'*^" < 
< [w]f + Hf [u;t]f < (1 + [i*])N^^'"- 

Proof Let u,t > and w £ L"(i?). On the one hand, for x,x' E E and 
s < t*{x) A t*{x') At Au, one has 

|wt($(a;,s)) - wt($(a;',s))| = i - s), s)) - i - s), s)) . 

We now use the fact that w G L"(i?) which yields since s < u 

\wtmx, s)) - WtiHx', s))\ < Hf - s) - {x',t- s)\ = Hf " ^1- 

Hence, [t«t]f < [wjf and similarly one obtains [wt]^'*^" < [wjf + [w]^'". 
On the other hand, for x,x' G E such that t*{x) V t*{x') < t Au, one has 

\wtmx,t*{x)))~wtmx',t*{x')))\ 

= w(^^{{x,t-t*{x)),t*{x)))-w(^^{ix',t-t*{x')),t*{x'))) 
= w(^^{{x,t-t*{x)),t*{x,t- t*{x)))^ 
-w(^^{ix',t-t*{x')),t*{x',t-t*{x')))^ 

moreover since w E L"(_E) and t*{x,t — t*{x)) V t*{x',t — t*{x')) < u one has 
\wmx,t*{x)))~wmx',f{x')))\ < [w]f^^\ix,t-t*{x))~{x',t-t*{x'))\ 



18 



We conclude thanks to the Lipschitz continuity assumption |2.10| on t* provid- 
ing \{x,t — t*{x)) — (x',t ~ t*{x'))\ < (l + — x'\. Eventually, one has 

l_E,tAti ^ ^^..^E,u|■^ , rj-*l\ „„j ... ^ rt/\u/-n^\ |— | 



K]f < N* '"(1 + [t*]) and wt e V/^^iE). 



The following projjosition proves that Assumption |2.13| holds for the time 
augmented process {X)t>Q. 

Proposition 5.3 For all w G L"(£'), one has 

1. for alls,, S,' eE andse [0, ? (^) A Ait], 

w{^{S,s))^Qw{^{e,s)) \ < ([g] vi)Nf"|e-r|, 

2. for all t C &E such that V?(C') < u, 

Qw{^{s,t*{o))^Qy^{^{eJ*{a)) 

< ([Q] V 1) (1 + r ]) (Hf " + Mf") - 



in other words, Assumption 2.13 is satisfied with [Q] — ([Q] V + [i*]). 



Proof As in the previous lemma, for all t > 0, we will denote Wt the function 
of B{E) defined by wt — w{-, t). For ^ — {x,t) £ E and w £ L"(£'), one has, by 
definition of Q, 

QwiO^ [ ^w{OQ{{x,t),de)= f w{z,t)Q{x,dz) ^Qwt{x). (6) 

We may now check the regularity assumption on Q. Let £, = {x,t) and 
S' = {x',t') e E. Let s € [0;?(C) A ?(^') A u]. Thanks to the definition of $ 
and equation (|6| one has 

)w{^{i,s))-Qw{^{i',s)) 

= s), t + s) - s),t' + s) 

= \Qwt+s{Hx,s)) - Qwt'+s{<^{x' ,s))\ 

We split it into the sum of two differences 

\Qwt+s{<^{x,s)) - Qwv+s{Hx' ,s))\ 

< |Qwf+,($(x,s)) - Qwt+s{Hx',s)) \ + \Q{wt+s - wt'+s){Hx' ,s))\ 



5.2 



wt+s e L 



(t+s)A« 



(E) so 



On the one hand, we recall that thanks to lemma 
that, since s < (t + s) Am, we may use the Lipschitz continuity assumption 2.13 
on Q and the first term is bounded as follows 



\Qwt+s{'^ix,s)) -Qwt+si'^ix^s))} < [Q][wt+s]i 



E,(t+s)Au, 
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5.2 



also provides [wt 



-,E,{t+s)Au 



+s\l 



< On the other hand, 



Moreover, lemma 
and more basically, the second term in the above equation satisfies 

\Q{wt+s-Wt'+s)mx',s)) \ < [w]f'y~t'\. 
Eventually, one has 

)w{^{^,s))-Qw{^{e,s)) < ([Q] Vl)Hf'1?-r| 



We now reason similarly to bound 



where ^ = {x,t) and C = {x',t') e E are such that t*{^) V t*(^') < u. Equa- 
tion ^ yields 

= \Qwt+t'ix){'^{x,t*{x))) -Qwe+fix'){^ix'^t*ix')))\ 

that we now spilt into 

I Qwt+t' (x) {<i>{x, t* (x))) - Qwf+t' (x') {<^{x', t* {x'))) I 
< \Qwt+t'{x){'^{x,t* {x))) - Qwt+t'[x){^{x\t* {x')))\ 

+ \Q{Wt+t'(x) - QWt'+t'(x')){'^{x' ,t*{x')))\. 



Thanks to lemma 



5.2 



Wt+f(x) € Lc*^* (^))^"^^^^ Moreover, we assume, with- 
out loss of generality that t*{x) > t*{x') so that t*{x) \/t*{x') < {t + t*{x)) Am. 
Therefore, the first term in the a bove equation is boun ded, thanks to the Lips- 

by [Q]((l + r])Nf" + 



chitz continuity assumption 



2.13 



on Q and lemma 



5.2 



[w]f '")|a; — x'\. More basically, the second term is bounded by [w]f '"|t — t' + 
t*{x) ~t*{x')\ < Hf'"(|i-i'| + [t*]\x~x'\). Eventually, one has 

Qw;($(C,?(e)))-Qu;($(e',r(C'))) 

<[g](i + r])Nf "la; - x'\ + Mf - A + \t- 1'\ + [t*]\x - x'\) 
< ([g] V 1) (1 + ri) (Nf « + Nf ") - n. 

Hence the result. □ 

Consequently, we may apply our numerical method to the time augmented 
process {Xt)t>o- In other words, for / g LdE), c £ Lc{dE) and £ E, our 
approximation scheme may be used to compute 



N 



(7) 



We will now see that the time depending functional and the deterministic time 
horizon problems boil down to computing such quantities JAr(/, c)(^) for suitably 
chosen functions I and c. 
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5.2 Lipschitz continuous cost functions 



We first consider the time depending functional problem with Lipschitz contin- 
uous cost functions. Thus, let then I S Lc(£'), c G hddE) and x e -E, we wish 
to compute 



E, 



/ l{Xt,t)dt + y^c{X^-,Tj)l^x ^edE} 
Jo f-{ 



It is straightforward to notice that this quantity may be expressed using the 
time augmented process starting from the point = {x, 0). Indeed, one has 



/ l{Xt,t)dt + y2c{X^-,T,)l{x ^edE} 
Jo _i 



where Jn{1,c){^o) is given by equation ([t]). Although they are time depending, 
the cost functions I and c are seen, in the left-hand side term, as time invariant 
functions of the time augmented process. The expectation of the time depend- 
ing functional is therefore obtained by computing the expectation of a time 
invariant functional for the transformed PDMP thanks to the approximation 
scheme described in Section |4] This is what expresses the following theorem, 
which proof stems from the previous discussion. 

Theorem 5.4 Let I G Lc(-EO '^^d c G "LddE) and apply the approximation 
scheme described in SectionVAto the time augmented process {Xt)t>Oj one has 
then 



JV 



l{Xt,t)dt + Y.''^^T7^Tj)^ 



Vo 



< eN{l,c,X,A). 



where we denoted SNil, c, X , A) the bound of the approximation error provided 
by Theorem \4.5\ when our approximation^ scheme is applied with cost functions 
I and c to the time augmented process {Xt)t>o. 

Remark 5.5 The quantity 6^(1, c, X, A) is compute d wi th respect to the process 
{Xt)t>o instead of {Xt)t>o as presented in Theorem so that one has 

e^{l,c,X,A)= J2 (2[vk+if\\Zk+i-Zk+i\\p 



k=0 



+ ([F]'^ + A[F]'^)\\Sk+i~^k+i\\p 



A 



where (Zk, Sk)keN denotes the sequence of the post-jump locations and the inter- 
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jump times of the time augmented process (Xt)t>a and with 



1^ ' , 



[F][ - Cfim 
{F]i = [tie,, 
[F]2 = Cu 
[F]'i = a, 

K]f <eC'**CA (^K{A,v,,^^) + nCt.[\]i(Ct>Ci + Cc))+Ct.[l]t 

K]f < e^'*^- {Ct.CiCx + 2Ci + CxC + {2n ~ l)Cx{Ct.Ci + C^) + Cu 

Mf <K]f + r]K]f. 

and for all w G l^c{E), K{A,w) = Ei + E2A + E^[w]f + E^Cyj + [Q][w\'^ where 
eventually 

[Q]= {[Q]yi){i + [t*]), 

El = 2[l]fCt. + Ci{[t*] + 2CUM1) + [c]f (1 + Ct.Cx) 

+Cc{2[X]iCt> + CxCl[X]i + 2[t*]Cx), 
E2= C,Ct'Cx[t% ^ 
E^= {\ + Ct'Cx)[Ql 
E^^ 2Cx[t*] + Ct'[X]i{2 + Ct>Cx). 

5.3 Deterministic time horizon 

In the context of applications, it seems relevant to consider a deterministic 
time horizon tf. For instance, one may want to estimate a mean cost over a 
given period no matter how many jumps occur during this period. Actually, 
we will choose a time horizon of the form tf l\ Tn with TV large enough to 
ensure the A''-th jump will occur after time tj with a high probability i.e. that 
Pa; (TV < tf) be close to zero. For a discussion concerning the choice of such 
N, and in particular a theoretical bound of the probability F^iT^ < tf^, we 
refer to [3]. Simply notice that in practice, this probability may be estimated 
through Monte-Carlo simulations. We thus intend to approximate the following 
quantity for I g Lc{E), c € LddE) and x € E: 



Tjv At/ 



Tjv 



l{Xt,t)dt+ J2 <XT7,Tj)l{x^_eOE} 
N 

l{Xt,t)l{t<tf}dt + ^c{Xrj.-,Tj)l{x^-edE}'^{T,<ti 



The natural approach would consist in killing the process at time tf as M.H.A. 
Davis suggests in j5|, section 31, and applying our method to the new process. 
However, the killed process will not necessarily fulfill our Lipschitz continuity 
assumptions because of the discontinuity introduced at time t f . 
A second idea would then be touse the previous results, to consider the time aug- 
mented process, and to define ^(x, t) — ^(x, t)l{t<tj,} andc(a;,i) — c(a;, i)l{t<t^j.. 
However, a similar problem appears. Indeed, such functions I and c are not Lip- 
schitz continuous and our numerical method requires this assumption. In the 
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rest of this section, we will see how to overcome this drawback. On the one 
hand, we prove that the Lipschitz continuity condition on I may be relaxed so 
that our numerical method may be used directly to approximate Jn(1,c) for 
any c G hddE). On the other hand, in the general case, we will deal with the 
non Lipschitz continuity of c by bounding it between two Lipschitz continuous 
functions. 



5.3.1 Direct estimation of the running cost term 

Let us explain how the Lipschitz continuity condition on the running cost func- 
tion may be relaxed so that Theorem |4.5[ stating the convergence of j3ur ap- 
proximation scheme, remains true when the running cost function is l{x,t) — 
l(x, i)l{(<t^.} with I e Lc(-E) and the boundary jump cost function is c € Lc(9i?) 
(although with slightly different constants in the bound of the convergence rate) . 
Indeed, the running cost function I appears inside an integral that will have a 
regularizing effect allowing us to derive the required Lipschitz property of the 
functional in spite of the discontinuity of I. Details are provided in Appendix [b| 

Consequently, our approximation scheme may be used to compute Jnih c)(C) 
for any ^ € E. We recall that Jjv is defined by ([t]) and that for all x £ E, one 
has 

Jn(1,c){x,0)=-E: 

We now turn to the indicator function l{7i^.<j^.} required within the boundary 
jump cost term. 



l{Xt,t)dt + y2c{X^-,T,)l{x _6aB} 

3 = 1 



5.3.2 Bounds of the boundary jump cost term 

We explained how the Lipschitz continuity condition on I may be relaxed. How- 
ever, when it comes to c, this condition cannot be avoided and our numerical 
method cannot be used directly with c{x,t) = c{x,t)l^t<tf}- We overcome this 
drawback by using Lipschitz continuous approximations of the indicator func- 
tion. Indeed, for S > 0, we introduce the real-valued functions and ub 
defined on M by 

= < -B{t - tf) if tf - i <t< tf, 




iit<tf, 

UB{t)={ -B{t~tf) + \ iitf<t<tf + ^, 



if tf + j^< t. 



The following lemma is straightforward. 



Lemma 5.6 For all t > 0, lim^^+oo ^bW — l[0;t/[(^) "^^^ limB-j.+oo ub(^) = 
l[0;t/](^)- Furthermore, for all B > 0, Ub CLuduB are Lipschitz continuous with 
Lipschitz constant B. Eventually, one has \ub ~ -"-[O:*/]! — 1; ~ l[0;t/]| ^ 1 
and 

Ub < l[0;t/] < UB- 
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Thus, define for I e LdE) 

T{x,t) = l{x,t)l^t<tj} (8) 
and for c € Lc{dE) and for all B > 0, 

Cg{x,t) = c{x,t)u^{t) and CB{x,t) = c{x,t)uB{t). (9) 
We now check that these functions satisfy our Lipschitz continuity conditions. 

Proposition 5.7 The functions Cg andcs belong to LcidE) w'tft [c^]*, [cs]* < 
[cl + BCc{lV[t*]). 

Proof We prove the result for c^, the other case being similar. For all ^ = 

{x,t),^' = {x',t') e E, one has 

\cB{H^,t*{i)))-CB{^S,',t*{a))\ 

= H^ij*m^B{t+t*{0)-<^i\t*{a))uB{t'+r 

< [c]*|e - ^'1 + Ce|MB(* + t*m - Ms(i' + *~*(^'))| 

< [c]*|e - <e'l + CcB{\t - t'\ + [t*]\x - x'\) 

< {[c]^ + C,B{iy[e]))\i-^'\. 

Hence the result. □ 

Therefore, the functions and cb are acceptable boundary jump cost func- 
tions and we may bound the deterministic horizon expectation by 



Jn{1,Cb){x,^) < E, 



Jo ^ 3 



< Jn{Icb){x,0). 



The following proposition provides the convergence of the bounds. 
Proposition 5.8 For all x G E, one has 

lim Jn{1,Cb){x,0) = lim Jn{1,cb){x,0) 



/ l{Xt,t)dt + Y c(^Tr > Ti)l{x^_ eaB} l{T,<t,} 



The above convergence holds for every tf > in the case of Jn{1,cb){x,0) but 
only for almost every tf > with respect to the Lebesgue measure on R in the 
case of Jn{1,Cb){x,0). 
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Proof Let x G E. We first consider Jjv(Z, cb)(x, 0). 



N N 




TV 



< ^1 c{X^- ,Tj)\\l <t, } - UB (T, ) 




where tp^ is the distribution function of Tj. For all j < N, limB-s-+oo + 



j^) — LPj{tf)j — since (^^ is right-continuous which shows the required conver- 
gence. 

We now turn to the case of JN{l,Cg){x, 0). Similar computations yields 



One cannot conclude as in the previous case since cpj is not necessary left- 
continuous. We therefore assume that tf is not an atom of any of the laws of the 

random variables Tj. Then, for all j < N, limB->.+oo (^fjitf) ~ ^j{if ~ = 
and the result follows. Indeed, the set of the atoms of Tj is at most countable 
so that the convergence holds for almost every tf w.r.t. the Lebesgue measure 
on M. □ 

5.3.3 Bounds in the general case 

The previous results show that the deterministic horizon expectation may be 
bounded by applying our numerical method with I and successively Cg and cb- 
In other words, we have shown the following theorem: 

Theorem 5.9 Let I e Lc(-E) and c E hc{dE). Let B)o<k<N (respectively 
(^fe.s)o<fc<JVy' be the sequence of random variables (yk)o<k<N described in Sec- 
tion^when applying our approximation scheme to the time augmented process 
(Xt)t>o with cost functions I and (resp. cb) defined by Q and The 
bounds of the approximation error provided by Theorem\4.5\ are respectively de- 



N N 




N 



<Ccj2{iPjitf)~ip,itf--)y 
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jted 8^(1, Cq, X , A, B) and £n{1,cb,X,A^B). One has the 



cTNMf 



<E, 



Jo _i 



< Vo^b + sn{Icb,X,A, B). 

Remark 5.10 In the previous theorem, the quantity 6^(1, Qb^ X , A, B) (and 
similarly e]\[{l,CB, X, A, B)) is computed with respect to the process {Xt)t>a in- 
stead of {Xt)t>o as presented in Theorem 4-5 so that one has 



N-l 



eN{l,CB,X,A,B) = h[vk+if\\Zk+i - Zk+i 



fc=0 



{2[v,f + [F][ + [f]';a+[f]';'b)\\z, ~z,\\, 

+ {[F]'^ + [F]'^A)\\Sk+i-S,+i\\p 



A 



where (Zk, Sk)ki£N denotes the sequence of the post-jump locations and the inter- 
jump times of the time augmented process {Xt)t>o and with 



[F]r = Ce(ivr]), 

[Vn]f < e^-'^- (k{A, B, Vn-l) + nCt> [\]l(Ct^Cl + Ce)) + Ct> \ 



and for allw G L^iE), K{A,B,w) = E[+E'^ B+E2A+E3[w]f +E4C^+[Q][w]l 
where eventually 

E[ = 2Wf a. + Ci{[t*] + 2C2.[A]i) + [c]f (1 + ^.Ca) 

+Cc(2[A]iCt. + CxCl[X]i + 2[t*]Cx), 
E'{= C,{ly[t*]){l + Ct-.Cx) 



The other constants remain unchanged and we refer to remark |5.<5| for their 
precise expressions. 

Furthermore, it is important to stress the fact that applying twice our nu- 
merical method does not increase significantly the computing time. Indeed, the 
computation of the quantization grids is, by far, the most costly step. These 
grids, that only depend on the dynamics of the process, may be stored off-line 
and used for the approximation of both bounds. 



The choice of B. We now discuss the choice of the parameter B, the dis- 



cussion is quite similar to the one concerning the choice of A in Section 4.2 



proposition |5.8| suggests that B should be chosen as large as possible. How- 
ever, choosing a large value for B will lead to large Lipschitz constants that 
will decrease the sharpness of the bounds £Ar(^,c^,X) and £iq{l,CB, X) for the 
approximation error provided by Theorem |4.5[ Indeed, it is easy to check that 
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[vn] grows linearly with B (see the precise expressions of the Lipschitz constants 
above). Thus, in order to control the error proposed by Theorem 4.5 
essary that the order of magnitude of the quantization error \\Q. 
most 4. 



e 



it is nec- 
JL be at 



6 Numerical results 

6.1 A repair workshop model 

We now present a repair workshop model adapted from [S], section 21. 

In a factory, a machine produces goods which daily value is r(x) where 
X € [0; 1] represents a parameter of evolution of the machine, a setting chosen 
by the operator. For instance, x may be some load or some pace imposed on the 
machine. This machine, initially working, may break down with age-dependent 
hazard rate A(t) and is then sent to the workshop for reparation. Besides, the 
direction of the factory has decided that, whenever the machine has worked 
for a whole year without requiring reparation, it is sent to the workshop for 
maintenance. The daily cost of such a maintenance is q{x) while the daily cost 
of a reparation is p{x), with reasonably p{x) > q{x). We assume that after a 
reparation or a maintenance, that both last a fixed time s, the machine is totally 
repaired and is not worn down. 

We therefore consider three modes: the machine is working (m — 1), being 
repaired (m = 2), undergoing maintenance (m = 3). The state of the process at 
time t will be denoted Xt — {mt, Ct: t) where Ct is the time since the last change of 
mode (this component is required since the hazard rate A is age-dependent) . The 
state space isE= {{l}x [0; 365] x M+) U ({2} x [0; s] x M+) U ({3} x [0; s] x M+) . 
In each mode, the flow is <I>m((C,i),'") — {( + u,t + u) . Concerning the transi- 
tion kernel, from the previous discussion, one may notice for instance that from 
the point (1, (, t), the process can jump to the point (2, 0, t) \i C, < 365 and the 
jump is forced to (3,0,t) if C = 365. Figure [l] presents the state space and an 
example of trajectory of the process. 

Our aim is to find the value of the setting x that maximizes the expected 
total benefits B{x) i.e. the discounted value (the interest rate is denoted p) of 
production minus the costs related to maintenance and reparation over a period 
tf = 5 years: 



where 



B* = sup B{x) 



"(a;)l{,„,=i} -p{x)l{,nt=2} - q{x)l{rnt=3})dt 



We will use the following values r(x) = x, p{x) = 100a;2, q{x) = 5, s 

0.03 

365 



7 days, p = and A represents a WeibuU distribution with parameters a — 2 



et 13 = 600. 
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Figure 1: An example of trajectory. The process starts from the point Zq in 
mode m = 1 (machine in service). The machine may be sent to the workshop 
for repairs (m=2) or maintenance (m=3). 



Our assumptions clearly hold so that we may run our numerical method. 
We first need to find A' G N such that P(i,o,o)(TAr < */) be small. Monte-Carlo 
simulations lead to the value iV = 18. For a fixed x E [0; 1], we will therefore 
compute Jn{1,0){1,0,0) where l{m,C,t) = e~P* {r{x)l{rn=i} - p(a;)l{m=2} - 
q{x)l^rn=3})^{t<tf}- Eventually, notice that we could have chosen r, p and q 
slightly more generally by allowing them to be time-dependent. 

It is important to stress the fact that, once the Markov chain associated 
to the process is quantized, we will be able to compute the approximation of 
B{x) almost instantly for any x G [0; 1] because the same grids are used for 
every computation. Thanks to this flexibility, we are able to draw the function 
X — >■ B(x) and, thus, to solve the above optimization problem very easily. This 
is a very important advantage of our method. Indeed, if we computed B{x) 
through standard methods such as Monte Carlo simulations, we would have to 
repeat the whole algorithm again and again for each value of x and solving the 
optimization problem would be intractable. 

The following figure represents the approximation of the function B com- 
puted on a constant step grid of [0; 1] with step 10~^. This leads to the solution 
of the earlier optimization problem. Indeed, we obtain B* = B{x*) = 537.84 
where x* = 0.78 is the value of the setting x that maximizes the benefits of the 
factory. 

Let now x = 0.78, the following table presents the values of Vjv, which are the 
approximations of B{x), for different number of points in the quantization grids. 
A reference value is obtained thanks to Monte-Carlo method (10^ simulations) 

BMonte-Carlo = 537.69. 
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Figure 2: The function B drawn with 500 points in the quantization grids. 



Points in the quantization grids 


Vn 


relative error to 537.69 


20 points 


542. U 


0.83 % 


50 points 


539.57 


0.35 % 


100 points 


538.24 


0.10 % 


500 points 


537.84 


0.03 % 



From a computational time point of view, we already explained that the 
computation of large quantization grids is, by far, the most costly step since 
it may take up to several hours whereas the approximation of the expectation 
that follows is then almost instantaneous. However, we may notice, in the 
above table, that grids containing only 50 points yield a quite accurate result 
with merely 0.35 % error. Such grids only require a few minutes to be designed. 



Remark 6.1 We already noticed that the same grids may serve several pur- 
poses. For instance, we may also have been interested in the computation of 
the mean time spent by the machine in the workshop by taking l(m,(^,t) ~ 
1 r 1 • 

|me{2;3}| 

6.2 A corrosion model 

We consider here a corrosion model of an aluminum metallic structure. This 
example was provided by Astrium. It concerns a small structure within a strate- 
gic ballistic missile. The missile is stored successively in three different envi- 
ronments which are more or less corrosive. It is made to have potentially large 
storage durations. The requirement for security is very strong. The mechanical 
stress exerted on the structure depends in part on its thickness. A loss of thick- 
ness will cause an over-constraint and therefore increase a risk of rupture. It is 
thus crucial to study the evolution of the thickness of the structure over time. 
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Let us describe more precisely the usage profile of the missile. It is stored 
successively in three different environments, the workshop (denoted to = 1), the 
submarine in operation (to = 2) and the submarine in dry-dock (to — 3). This is 
because the structure must be equipped and used in a given order. Then it goes 
back to the workshop and so on. The missile stays in each environment during 
a random duration with exponential distribution. Its parameter Am depends 
on the environment. The degradation law for the thickness loss then depends 
on the environment through two parameters, a deterministic transition period 
77m and a random corrosion rate p uniformly distributed within a given range. 
Typically, the workshop and dry-dock are the most corrosive environments but 
the time spent in operation is more important. The randomness of the corrosion 
rate accounts for small variations and uncertainties in the corrosiveness of each 
environment. 

In each environment to G {1; 2; 3}, the evolution over time of the thickness 
loss dm satisfies: 

dm{p,s) = p{s + ^ra[er^ ■ (10) 
Table [T] gives the values of the different parameters. 



environment 1 environment 2 environment 3 



A™ (h-^) 


(17520)-! 


(131400)-! 


(8760)-! 


Tim (h) 


30000 


200000 


40000 


p (mm.h"!) 


[10-6,10-5] 


[10-^ 10-6] 


[10-6,10-5] 



Table 1: Numerical values of the parameters of the corrosion model. 

Initially, the structure is in environment to = 1 and the thickness loss is null. 
One draws the corrosion rate po uniformly distributed in the interval [10-6, 10-5] 
and the time of the first change of environment T\ exponentially distributed 
with parameter Ai = (17520)-! hours-!. The corrosion starts according to 



Equation (10 1 so that, for all < t < Ti, the loss of thickness is di(pQ,t). The 
structure then moves to environment 2 and the process restarts similarly: a new 
corrosion rate is drawn according to an uniform law on [10-^, 10-6], the time 
of the second jump T2 is drawn so that T2 — Ti is exponentially distributed with 
parameter A2 — (131400)-! hours-! and for Ti < t < T2, the loss of thickness 
is c?i(po, Ti) + d2{pTi ,t — Ti) and so on. 

At each change of environment, a new corrosion rate p is drawn according 
to a uniform law on the corresponding interval. The thickness loss, however, 
evolves continuously. 

We are interested in computing the mean loss of thickness in environment 2 
until a given time i/ = 18 years. 



Modelization by PDMP 

The state space E. The loss of thickness will be modelized by a PDMP 
whose modes are the different environments. Let then M ~ {1,2,3}. The 
PDMP {Xt)t>o will contain the following components: the mode to G M, the 
loss of thickness d, the time since the last jump s (this is to ensure that the 
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Markov property is satisfied), the corrosion rate p and the time t (since we 
consider the time-augmented process) . Notice that clearly, one has always s <t 
so that we reasonably consider the following state space: 

E= {{m, d, s,p,t) e M xR+ xR+ X [10"'^; 10"^] x K+ such that s < t} . 



The flow The flow is given for all u > by 





m 


\ 




/ 


m 




d 








d + dmip,s + u) 




s 




,u) = 




s + u 




P 








p 


V 


t 


I 




\ 


t + u 



I 



The transition kernel Q. Let us now study the jumps of this process. When 

the process jumps from a point x = (m, d, s, p, t) G E, m becomes m+l modulo 
3 (denoted m+l [3]), d and t remain unchanged, s becomes 0. Only p is randomly 
drawn, according to a uniform law on an interval [pmin\ Pmax] that depends on 
the new mode. One has then for w G B{E), x = (m, d, s, p, t) G E, and u>0, 



( m\ 

d 
s 
P 

\ t J 



i)) = Qi 



f m \ 

d + dm{p, s + u) - dm{p, s) 
s + u 



\ 



p 

t + u 



J 



-'max ymzn 



W 



( m+l[3] \ 

d + dm{p-, s + u)- dm{p, s) 




V 



p 

t + u 



dp. (11) 



The cost function I. The function I G B{E) will be the cost function to 
compute the mean loss of thickness in mode 2. It is defined as follows: for all 
X = (m, d, s, p,t) G E and u>0 



imx, u)) = p(l-^e ^ 1{„=2} = ^ {dm{p, S + U)) l{m=2}- 

One then defines l{^{x,u)) = l{^{x,u))l^t^u<t^} so that, one has 

L{x,u) = f l{^{x,u'))du' 
Jo 

P-t)+ 

l{<^{x,u'))du' 

(dm{P, S + UA{tf - t) + ) - dm{p, S)) l{m=2} 



(12) 



= 1 



that is indeed the thickness lost in mode 2 from the point x = {m,d,s,p,t) 
during a time uA{tf—t)~^. 
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The assumptions. Assumption 2.1 and 2.9 are clearly satisfied. Moreover, 
it is straightforward, from (12 1, to check that I £ Lc(-E) so that Assumption 3.1 
holds. 



We now turn to Assumption |2.13| and we will see that, although it does 
not hold for any function w G L"(i?), it holds for a sufficiently big sub class of 
functions. We first need to make a remark. Recall that for all x — (to, d, s, p, t) e 

E and for all k G {0, N}, one has UAr_j.(x) = Ej. J^'° l(^{x, u)) l{t+„<j^.}0?u . 

Therefore, for all k € {0, ■■■,N} the function Vk as well as the function / satisfy 
the following condition: 

for all X = (to, d, s, p,t) £ E such that t > tf, one has w{x) — 0. (13) 



The next step consists in proving that Assumption 2.13| although it is not 



satisfied for any function w G ^^{E), holds for any function w £ J-^ciE) that also 



satisfies condition (13). This is done in Lemma 6.2 and it is sufficient because 



in the proof of the theorem that ensures the convergence of our approximation 



scheme. Assumption 2.13 is only used with the functions {vk)ke{o,...,N} that do 
satisfy condition (13 1. 



Lemma 6.2 There exists [Q] E M+ such that for all v >0 and w G Lc('^) ^^^^ 
satisfies condition (|13[), one has for all x, x' (£ E and < u < v, 



\Qw{<P{x,u)) ~ Qw{<^>{x' ,u)) \ < [Q][w]f'''\x-x'\. 

Proof Let x = {m,d,s,p,t) and x' = {m' ,d' ,s' , p' ,t') G E with for instance 
t < t' . First we may choose m — m' , otherwise, jx — x'| = +oo and there is 
nothing to prove. Now, we are facing three different cases: 

9 a t f < t + u < t' + u, then one has Qw{^{x,u)) = Qw{^{x' ,u)) = 
because w satisfies condition (|13[) and there is nothing to prove. 



if t+u <tf< t'+u, notice that Qw[^{x', u)) = Qw(^{{m', d' , s' , p' ,tf),u)) 
(this stems from condition ( [l3| ) so that we are reduced to the following 

case, 

We assume from now on that t + u < t' + u < tf. We now intend to 



bound \Qw{^{x,u)) — Qw{^(x\u))\. It is clear from equation (111 that 



we only need to prove that the function (p, s) — > dm{p, s), defined by ( 10 1 
is Lipschitz continuous w.r.t. both its variables on the set [10~^; 10~^x 
[0; tf]. Indeed, we have s < t and s' < t' so that s, s', s + m, s' + m < tf. 
Standard computations yield: 

3 

\dm{P,s) - dm{p' ,s')\ < S\p- p'\+ -p'\s ~ s'\ 

< tf\p-p'\ + ^10-'\s~s'\. 
Hence the result. □ 



Eventually, Assumption |2 . 1 1 is not satisfied because in our corrosion model, 
one has t*{x) = +oo for all x G E. Besides, we may notice that the previous 
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proof would have been more straightforward if t* had been bounded. Indeed in 
that case, we would have had s, s' , s + u, s' 



condition ( 13 ) would have been unnecessary. 



u < Cf and the introduction of 
Nevertheless, we have been able 
to overcome the drawback of having t* non-bounded by noticing that somehow 
the deterministic time horizon tf plays the part of the missing Cf. This is the 



meaning of condition (13): roughly speaking, we do not consider what happens 
beyond tf. 

More generally, we will now see that in our deterministic time horizon problem, 
the boundedness of t* may be dropped and our re sults remain true replacing 
Cf by tf. This is clear in the case of Proposition A. 2 because the function I 
satisfies the condition (13 1. Proposition A. 7 remains also true replacing Cf by 
tf. Indeed, on the one hand, it is clear that L(x, u) < tfCi. On the other hand, 
when computing |w„(<i>(a:, u)) — w„(<i>(x', we are facing three different cases 
(as in the proof of Lemma 6.2): 



a tf < u < u', one has w„(<i>(x, m)) 
condition ( 13 1), 



Vn{^{x' ,u')) — (this stems from 



• if u < tf < u' , one has |w„(<i>(x, m)) — Vn{^{x' ,u'))\ ~ \vn{^{x,u)) — 
Vn{^{x' ,tf))\ since u„(<i>(a;', u')) = Vn{^{x' ,tf)) — (condition (13) once 
again) so that we are reduced to the following case, 

• If u < u' < tf, the computations remain unchanged and tf replaces Cf 
as a bound for u and u' . 



Numerical results 

Table [2] presents the values of the loss of thickness in environment 2 obtained 
through our approximation scheme as well as a Monte Carlo approximation 
(obtained with 10^ simulations) and the relative errors of our values w.r.t. the 
Monte Carlo value. Figure |3] presents respectively the empirical convergence 
rate. The convergence rate, estimated through a regression model is -0.35. 
This is roughly the same order of magnitude as the rate of convergence of the 
optimal quantizer (see for instance [9 ) since here the dimension is 3 (indeed, m 
is deterministic and s = immediately after a jump so that we only quantize 
the variables p, d and t). 

Besides, Table[3]presents the CPU time to compute the expectations from the 
quantization grids (computations are run with Matlab R2010b on a MacBook 
Pro 2.66 GHz 17 processor). It can be seen that, once the quantization grids are 
obtained, our approximation scheme is performed very quickly. 



7 Conclusion 

We have presented an efficient and easy to implement numerical method to 
approximate expectations of functionals of piecewise-detcrministic Markov pro- 
cesses. We proved the convergence of our algorithm with bounds for the rate of 
convergence. 

Although our method concerns time invariant functionals, we proved that we 
are able to tackle time depending problems such as Lipschitz continuous time 
depending functionals or deterministic time horizon expectations. Indeed, we 
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Quantization grids 


^^0 


error 


20 points 


0.038386 


4.43 % 


50 points 


0.037804 


2.85 % 


100 points 


0.037525 


2.09 % 


200 points 


0.037421 


1.81 % 


500 points 


0.037264 


1.38 % 


1000 points 


0.037160 


1.10 % 


2000 points 


0.037041 


0.77 % 


4000 points 


0.037007 


0.69 % 


6000 points 


0.036973 


0.57 % 


8000 points 


0.036944 


0.49 % 


10000 points 


0.036911 


0.40 % 


12000 points 


0.036897 


0.36 % 


Monte Carlo 


0.036755 





Table 2: Approximation of the mean loss of thickness (in mm) in environment 2 
for different numbers of points in the quantization grids and a Monte Carlo 
approximation (10^ simulations). 




Figure 3: Logarithm of the error when approximating the loss of thickness in 
environment 2 w.r.t. the logarithm of the number of points in the quantization 
grids. The empirical convergence rate, estimated through a regression model, 
is -0.35. 



proved that, thanks to the introduction of the time augmented process, time 
depending problems may be seen, paradoxically, as special cases of the time 
invariant situation. 

Our method is easy to implement because it merely requires to be able to 
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Quantization grids 


CPU time (second) 


20 points 


0.0059 


50 points 


0.0085 


100 points 


0.014 


200 points 


0.034 


500 points 


0.12 


1000 points 


0.37 


2000 points 


1.5 


4000 points 


5.6 


6000 points 


13 


8000 points 


24 


10000 points 


35 


12000 points 


54 


Monte Carlo (10** simulations) 


W16000 



Table 3: CPU time. 



simulate the process. Furthermore, although the computation of the quanti- 
zation grids may be quite time-consuming, it may be performed preliminarily 
because the grids only depend on the dynamics of the process and not on the 
cost functions I and c. Therefore, they may be stored off-line and serve several 
purposes. As it is illustrated by the examples presented in Section [6] storing the 
grids provides to our approximation scheme efficiency and flexibility. Indeed, 
the computation of the expectation can be performed very quickly once the 
grids are available. Thus, if one decides for instance to modify the functional, 
the same grids may be used so that the new result is obtained very quickly. This 
flexibility is an important advantage over standard Monte-Carlo simulations. 
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A Lipschitz continuity of F, G and Vn 

The first lemma and the first proposition of this section present mainly the 
Lipschitz continuity of the functions and F. They are stated without proof 
because they are quite straightforward. 

Lemma A.l The function is Lipschitz continuous w.r.t. both its variables 
i.e. for all x, y Cz E and u, t £R, one has 

\S^{x,t)-S^{y,t)\<A[t*]\x-y\, 
\S^{x,t)-S^{x,u)\ <A\t~u\, 
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Moreover, one has for all x G E and t, s > such that t + s < t*(x), 

S^mx,s),t) = S^{x,t + s). 



Proposition A. 2 The function F, introduced in Definition \3.3\ is Lipschitz 
continuous w.r.t. both its variables. Forallx, y (z E andu,v € [0;t* (x) At* (y)], 
one has 



with 



\F{x,u) - F{y,v)\ < [Fl\x ~ y\ + [F]2\u - v\ 

[F]i =Ct4l]i + [c],+A[t*]C,, 
[Fh ^Ci + AC,. 



The two following lemmas are adapted from [6,, the second one being a 
special case of lemma A.l from Thus, they are stated without proof. 

Lemma A. 3 For h e Lc{E), {x,y) G E^ , and t < t*{x) At*{y) 



ft'ix) 



t'iv) 



h{<^{x,s))e-^^'''''Us- J h{<^{y,s))e-^^y''''>ds 

< {cAh], + {CUMi + [t*])CH)\x-y\. 
Lemma A. 4 For h G hc{dE) U Lc(i?) and x,y E, one has 

< (N* + Ch{Ct* [A]i + [t*]Cx)) \x-y\. 

We now introduce a definition that will be convenient in the sequel. For 
w e Lc{E), X £ E and t e [0;i*(a;)], we define 

Gtwix) = E, [{F{x, Si) + w{Zi)) l{s^>t}] 

= E, Si) + C{x, Si) + w{Zi)) l{s,>t}] ■ 

In particular, Gq = G. Since we know the law of (Zi, Si), it can be shown that 

Gtw{x) = Ti(x) + T2{x) + Tsix) + T^{x) + T^ix) (14) 



with 



Ti (x) = e-^(^'*) IqI°'^{x, s)ds, 
T^ix) /f ^""^ I o <i>{x, s)e^'^(^''*)ds. 

T3(x) = c o <^{x, t*{x)) jf'-''^ S'^ix, s)X o $(x, s)e-^(^''*)ds, 

T4(a;) = /f ^""^ (XQw) o <t>{x, s)e-^^^''Us, 

Tslx) = e-^(^'**(^»(Ow + c) o$(a;,r(a;)). 
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Proposition A. 5 For w G Lc(£'), {x,y) € and t e [0;t*{x) A t*{y)], one 
has 

\Gtw{x) - Gtw{y) \ < K{A,w)\x-y\, 
where K{A, w) ^ Ei + E2A + E3[w]i + E^Cw + [Q][w\^ with 

El = 2[/]iCt. + Ci{[t*] + 2Cl [A]i) + [c],(l + Ct^Cx) 

+C,{2[X]iCt*+CxCl[X]i+2[t*]Cx), 
E2^ C,Ct'Cx[t% 
Es^ {l + Ct'Cx)[Ql 
E^^ 2Cx[t*] + Ct4M'2 + Ct'Cx). 

Proof Let w G Lc(-E), {x,y) € E'^ and t e [0;t*{x) At*{y)]. From equa- 
tion (14l, we naturally split IGtW^x) — Gyw{y)\ into the sum of five differences. 

The first one is |Ti(a;) — Ti(j/)| and is bounded by 



|Ti(x) - Ti(2/)| < Ct.Ci e-^(-'*) - e'^^^^*^ 

< (cf,Ci[Xl + Ct,[l]i)\x-y\ 



J o ^{x,s) - I o ^{y,s)^ds 



The differences \T2{x) — T2{y) \ and |T4(x) — T4{y) \ may be bounded thanks 
to lemma A. 3 with successively h = I and h = XQw. Notice that Cxqw ^ GxC^ 
and [XQwjri Ca[Q]Ni + 



For the difference of the T5 terms, we use lemma |A.4| with h = Qw + c. 
Notice that Cqw+c < G^ + Gc and that [Qw + c]* < [Q] {[w]^ + ['w]i) + [c]*. 

Eventually, in order to bound 1X3(0;) — T3(y)|, we assume without loss of 
generality that t*{x) < t*{y) and we have 



\T,{x)~T,iy)\ 



< C, 



S^{x, s)X o <^{x, s)e"'^(^'") - (5^(y, s)X o s)e-^^y-'''^ 



ds 



t'iv) 



t*(x) 
t'(x) 



5^{y, s)X o s)e-^(y'^) ds + [c],Gt-.Gx\x - y\ 



I (^Gx\5^ix,s) -5^{y,s)\ + [X]i\x-y\ + Gx\e-^'^'=''^ - e-^(^^^)|) ds 
+G,[t*]Gx\x~y\ + [clCt.Gx\x~y\ 
< (CcGt* {CxA[t*] + [A]i + CAa.[A]i) + G4t*]Gx + [c].Gt^Gx) \x ~ y\. 

The result follows. □ 

The following lemma is stated without proof. It is indeed very close to 
Lemma 51.7 from [S]. 

Lemma A. 6 For all x E E and t E [0; t*{x)], one has 

«„($(x,t)) = e^("'*)Gtz;„_i(x)- / /o$(x,s)ds. 
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Proposition A. 7 For all n E {0, 1, A^}, f„ G Lc(£^) and one has 

a„ <n{Ct'Ci + Cc), 

Ml <e^"^^ (^K{A,Vn-i) + nCt4^]i(^Ct-.Ci + Cc))+Ct4l]i, 

Kh < e'^'*''^ (Ct^CiCx + 2Ci + C^C, + {2n - l)C^{Ct,Ci + C,)) + d, 

[Vn]* < [Vn]l + [t*][Vnh- 
K] < K{A,Vn-l), 



Proof Recall that for x E E, one has from definition |3. 3 1 

Vnix) = Gi;„_i(x) - E, Si)] + E, [Cix, Si)] + E,K_i(Zi)]. 
Thus, C„„ < CfCi + Cc + C„„_i < n{Ct'Ci + C^) by induction. 

Let us now turn to [vn]i. Lemma |A.6| yields 

\Vr^mx,t))~Vnmy,t))\ 

< \e''^^'''>GtVn^i{x)-e^^y^'^GtVn-i{y)\+ f \l o <i>{x,s) - I o ^y,s)\ds 

Jo 

< e^^-''^\GtVr,-i{x) - GtVn-iiy)\ + |Gtt;„_i(y)||e'^(^'*) - e'^^^'')! + Gt4l]i\, 

The result follows using proposition |A.5| and noticing that 

A{x,t)<Gt,C^, 
\Gtv^.i{y)\ < Gt.Ci + a + a„_, < n{Ct>Ci + Ce), 

We now turn to [w„]2- For x £ E and s, t £ [0,t*{x)] with s < t, one has 

|w„($(x,<)) -'y„($(a;,s))| 
<eA(-'*)|Gti;„_i(x)-G,^;„_i(:r)| + |G,^;„_i(a;)||e^(-'*)-e^(--^)| + Q|t-s|. 



Moreover, from (14), one has 



\GtVn^i{x) - GsVri-i(a:^)| 

< E, [\F{x,Si) + Vn-l{Zi)\l{,<s^<t}] 



< 



l{^{x,u))du 



,-A{x,s) 



l{^{x,u))du 



du 



+ \co^{x,t*{x))\ I (5^(a;,w)Ao$(x,u)e-'^(^'") 

du. 



du 



< 



(AQi;„_i)oci>(:E,w)e-^(--") 

(Gi.G,|e-^(-^*) - e-^(-'^)| + Ci\t - s\) + {Ci\t - s|) 
+ (c,Cx\t-s\^ + (GAC,„_Ji-s|). 
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and 



|gA(x,t) _gA(x,.)| <eC**CAC7A|i-s|. 



Eventally, the bound for [w„] is a direct consequence from proposition A. 5 □ 



B Relaxed assumption on the running cost func- 
tion 

In this section, we consider the approximation apphed to the time augme nted 
process so that the local characteristics a£e A and Q defined in Section 
Moreover, we consider a function I £ Lc(i?) and we define / G B{E) by 



5.1 



for all ^ = {x,t) e E, = l{x,t)l 



{*<«/}■ 



We intend t o pr ove that the convergence of^our approximation scheme, stated 
by Theorem 



4.5 



remains true if we choose I as the running cost function even 
though it does not fulfill the required Lipschitz conditions i.e. I ^ Lc(-E) . Indeed, 
the Lipschitz continuity of I is used four times in the proof of the theorem, once 
in proposition A.2[ twice in proposition A. 5 (when bounding the difference of 
the Ti terms and the one of the T2 ones) and once in proposition A. 7 (when 



bounding [vn]i)- In each case, the Lipschitz continuity of the running cost 
function I is used to bound a term of the form 



du 



for ^, C & E and s, s' G [0; t*{() A t*{^')], or of the form 



du 



(15) 



(16) 



for ^, ^' e and s £ [0; t*{£_) At*{C)] and where we naturally denoted A(^, u) — 
\{^{£_,v))dv. Concerning this second form, equation ( [T6| , notice that 



t*(«)At*(?') 



< 



du 



lo^(£^,u) -lo^(^',u) 



du 



< 



Ci 

At* («') 



du 



lo<i>{^,u) - lo u) du + CiCf, [A]i 

remains true with I as the running cost 



4.5 



so that, to ensure that Theorem 
function, it is sufficient to be able to bound terms of the form ( 15 ). This is done 
in the following lemma. 

Lemma B.l For £, ^ {x,t),^' = {x',t') e E and s e [0; ? (^) A F*(C')]; one has 
To^{^,u) -Jo^{C,u) du<{Ct4l]i + Ci)\^-ei 
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Proof Let ^ = {x,t),(,' = {x',t') E E a.nd s E [0;t*{£,) M* {(,')], one has 



Zo$(^,u) -;o$(^',u) 



< 



< 



/oci>(e,u)-Zoci>(^',,.) 



+ Ci 



l{t+ji<t/} ~ l{t'+«<t/} 



The left-hand side term is bounded by Ci* ~ CI since I e Lc(£'). For the 
right-hand side term, assume without loss of generality that t < t' , one has 



^{t+U<tf} ~ ^{t'+U<tf} 



l{t-t/<tl} — l{t'-t/<M} 



1 



{t-tf<U<t'~tf} 



so that the right-hand side term is bounded by Ci\t — t'\ < Ci\^ — C\- The result 
follows. □ 



Eventually, Theorem 4.5 remains true if we choose I as the running cost 
function. One only needs to slightly modify the Lipschitz constants given in 



propositions A. 2 A. 5 and A. 7 The terms Cf [l]i have to be replaced by Cf 



C Proof of Theorem 4.5 



The Lipschitz continuity of the functions Vk is proved by proposition |A.7[ Let 
now ^4 > and first notice that 

\Jn{1,c){x)- Vol < \Jn{1, c){x)- Vol + \Vo- Vol 



Proposition 3.2 states that \ Jn{1,c){x) - Vol < since Vo = J^{l,c){x). 



We now have to bound |Vb — Vo|. 

Some of the arguments of the proof are similar to the ones used in Theo- 
rem 5.1 from [B], thus we will not develop the details of the proof. Recall that 
II Vat — VnWp = and let fc G {0, — 1}. In order to bound the approximation 
error, let us split it into three terms || Vfc — Vfc||p < Si + S2 + S3 where 



^1 

^3 



\\vk{Zk) -VkiZk)\\p, 
\\Gvk+i{Zk) - Gk+iVk+iiZk)\\p, 
\\Gk+iVk+i{Zk) — Gk+iVk+i{Zk)\\p. 



The theorem is then a direct consequence from the three following lemmas, 
stated without proof, that provide bounds for each of these three terms. 

Lemma C.l The first term Si is bounded by 

\\vk{Zk) - VkiZk)\\p < [vk]\\Zk - Zk\\p. 
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Lemma C.2 The second term S2 is bounded by 

\\Gvk+i{Zk) - Gk+iVk+i{Zk)\\p 

< [vk+i]\\Zk+i - Zk+i\\p + {[vk] + [F]i)\\Zk - ZkWp + [F]2\\Sk+i - Sk+i\\p. 
Lemma C.3 The third term 'E.3 is bounded by 

\\Gk+iVk+i{Zk)-Gk+i%+i{Zk)\\p < bfe+i]||-^fe+i--^fe+i||p+||Vfe+i-Vfe+i||p. 
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